# Replication archive for:
# Galos and Coppock. 2023. 
# Gender Composition Predicts Gender Bias: A Meta-reanalysis of Hiring Discrimination Audit Experiments. 
# Science Advances. Forthcoming.

library(tidyverse)
library(metafor)
library(broom)
source("code/helpers.R")

survey_studies <- read_rds("data/survey_studies.rds")

survey_meta_fit <-
  rma.uni(
    yi = estimate,
    sei = std.error,
    method = "REML",
    data = survey_studies
  )


survey_studies$weight <- weights(survey_meta_fit)

survey_meta_ests <-
  tidy(survey_meta_fit, conf.int = TRUE) %>%
  mutate(
    study_id = "Meta-analysis",
    study_name_citation = "Meta-analysis",
    facet = "Meta",
    weight = mean(survey_studies$weight)
  )

# Join study level with meta analysis -----------------------

gg_survey <- 
  bind_rows(survey_studies, survey_meta_ests) %>% 
  ungroup() %>% 
  mutate(study_name_fct = fct_reorder(study_name_citation, estimate),
         entry = make_se_entry(estimate, std.error, 1),
         facet = factor(facet, levels = c("Study estimates", "Meta")))

# plot --------------------------------------------------------------------

g <- 
  ggplot(gg_survey, aes(estimate, study_name_fct)) +
  geom_point(aes(size = weight)) +
  geom_text(x = -22, aes(label = entry), hjust = 1, size = 2.5) +
  geom_linerange(aes(xmin = conf.low, xmax = conf.high), alpha = 0.3) +
  geom_vline(xintercept = 0, color = "red", linetype = "dashed") +
  geom_vline(data = select(survey_meta_ests, -facet), aes(xintercept = estimate)) +
  geom_rect(
    data = select(survey_meta_ests, -facet),
    aes(
      xmin = conf.low,
      xmax = conf.high,
      x = NULL,
      y = NULL
    ),
    ymin = 0,
    ymax = 100,
    alpha = .1
  ) +
  scale_size_area(max_size = 3) +
  theme_bw() +
  facet_grid(rows = vars(facet), scales = "free_y", space = "free") +
  coord_cartesian(xlim = c(-25, 25)) +
  theme(strip.background = element_blank(),
        strip.text = element_blank(),
        panel.grid = element_line(color = gray(0.95)),
        panel.grid.minor = element_blank(),
        axis.title.y = element_blank(),
        axis.text.y = element_text(size = 5),
        legend.position = "none") +
  labs(x = "Estimated average treatment effect in percentage points")

g